Valence Bond and von Neumann Entanglement Entropy in Heisenberg Ladders 
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We present a direct comparison of the recently-proposed valence bond entanglement entropy and 
the von Neumann entanglement entropy on spin 1/2 Heisenberg systems using quantum Monte 
Carlo and density-matrix renormalization group simulations. For one-dimensional chains we show 
that the valence bond entropy can be either less or greater than the von Neumann entropy, hence 
it cannot provide a bound on the latter. On ladder geometries, simulations with up to seven legs 
are sufficient to indicate that the von Neumann entropy in two dimensions obeys an area law, even 
though the valence bond entanglement entropy has a multiplicative logarithmic correction. 



Introduction.- Entanglement has arisen as a new 
paradigm for the study of correlations in condensed mat- 
ter systems. Measurements of entanglement between 
subregions, chiefly using entropic quantities, have an ad- 
vantage over traditional two-point correlation functions 
in that they encode the total amount of information 
shared between the subregions without the possibility of 
missing "hidden" correlations such as may occur in 
some exotic quantum groundstates. For example spin liq- 
uid states, where two-point correlation functions decay 
at large lengthscales, can exhibit topological order that 
is quantified by a "topological entanglement entropy" Q • 
This and other entropic measures are typically discussed 
in the context of the von Neumann entanglement entropy 
(S vN ), which for a system partitioned into two regions A 
and B, quantifies the entanglement between A and B as 



-Tr[pA In pa] 



(1) 



Here, the reduced density matrix pa — TrB\tp)(ip\ is ob- 
tained by tracing out the degrees of freedom of B. 

The properties of S' vN are well-studied in quantum in- 
formation theory. In interacting one-dimensional (ID) 
quantum systems, exact analytical results are known 
from conformal field theories (CFT); they show that, 
away from special critical points, S^ N scales as the size 
of the boundary between A and B. This so-called area 
law [H is also believed to hold in many groundstates 
of two dimensional (2D) interacting quantum systems, 
although exact results are scarce [J]. This has conse- 
quences in the rapidly-advancing field of computational 
quantum many-body theory, where it is known for exam- 
ple that groundstates of ID Hamiltonians satisfying an 
area law can be accurately represented by matrix prod- 
uct states 0. Tensor-network states and MERA give 
two promising new classes of numerical algorithms |6| 
that may allow simulations of 2D quantum systems not 
amenable to quantum Monte Carlo (QMC) due to the no- 
torious sign problem. However, these simulation frame- 
works are constructed to obey an area law; in order 
to be represented faithfully by them, a given 2D quan- 
tum groundstate must have entanglement properties also 



obeying the area law [j| • 

Unfortunately, entanglement is difficult to measure in 
2D, due to the fact that QMC does not have direct access 
to the groundstate wavefunction \ip) required in Eq. |T]). 
In response to this, several authors have introduced 
the concept of valence bond entanglement entropy (S ), 
defined for a spin system as 

^ B =ln(2)-Ad, (2) 

where Ma is the number of singlets (| jj,) — | J.T))/v2 
crossing the boundary between regions A and B. Unlike 
5 vN , S yB can be accessed easily in the valence-bond ba- 
sis projector QMC method recently proposed by Sand- 
vik [§]. As demonstrated in Refs. jail, S VB has many 
properties in common with S* vN , in particular the rela- 
tionship S\ B = S^ B , and the fact that S\ B = for 
regions "un-entangled" by valence bonds. A comparison 
of the scaling of S VB for (critical) ID spin 1/2 Heisen- 
berg chains shows good agreement with analytical results 
known from CFT, however in the 2D isotropic Heisenberg 
model it displays a multiplicative logarithmic correction 
to the area law 0, 0- If true also for S vN , this would 
have negative consequences for the simulation of the Neel 
groundstate using tensor- network simulations. 

In this paper, we compare S' VB calculated by valence- 
bond QMC to 5 vN accessible through density-matrix 
renormalization group (DMRG) simulations of the 
Heisenberg model on iV-leg ladders. For N = 1, the CFT 
central charge calculated by scaling 5 vN shows excellent 
agreement to c = 1, whereas S VB converges to c < 1. 
For N > 1, S* VB is systematically greater than S vN , a 
trend which grows rapidly with N. An examination of 
entanglement defined by bipartitioning multi-leg ladders 
shows a logarithmic correction for the valence-bond en- 
tanglement entropy, S* VB oc N In N (in agreement with 
Refs. [7, 8]), however for data up to N = 7, the von Neu- 
mann entanglement entropy convincingly shows a scaling 
of 5 vN oc N, thus obeying the area law. 

Model and Methods - We consider the spin 1 /2 Heisen- 
berg Hamiltonian, given by H = X!(ij> S» • Sj, where the 
sum is over nearest-neighbor lattice sites. Geometries 
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studied are ID chains, and multi-leg ladders with length 
L and number of legs N . We employ two complemen- 
tary numerical techniques, the valence-bond basis QMC 
and DMRG, both of which give unbiased approximations 
to the ground state of the Hamiltonian at zero temper- 
ature, is naturally accessible through the DMRG 
algorithm [lfj. At each step of the algorithm, the wave- 
function of the system is approximated by keeping only 
the states with largest coefficients in the Schmidt de- 
composition for a given bipartition into regions A and 
B = CA. To find the basis entering the Schmidt decom- 
position for region A, the reduced density matrix pA is 
calculated and diagonalized, thus allowing easy calcula- 
tion of Eq. ([1]). The truncation of the basis implies that 
only a lower bound for is calculated, so care must be 
taken to ensure that enough of the eigenvalue spectrum is 
included to converge S A N to sufficient accuracy; typically 
the number of states required is larger than that needed 
to converge the energy alone. We have kept up to 1800 
states per block for the largest ladders, with truncation 
errors under 10 -7 in all cases. 

S A B can be calculated using the valence-bond basis 
QMC proposed by Sandvik [||. The algorithm we use is 
the simple single-projector method, with lattice geome- 
tries constructed to match those used by DMRG. The 
ground state of the system is projected out by repeated 
application of a list of nearest-neighbor bond operators, 
a number of which are changed each step. The change is 
accepted with a probability of 2™ 6 /2" a where n a (rib) is 
the number of off-diagonal operators in the current (last 
accepted) step. Measured quantities such as energy or 
S\ B are then calculated by a weighted average over all 
the valence bond states obtained by this procedure. 

One- dimensional chain - We consider first the case of 
Heisenberg chains (N = 1) of length L, simulating both 
open (OBC) and periodic (PBC) boundary conditions. 
The DMRG algorithm requires the regions A and B to 
be topologically connected, so in ID the bipartition is 
defined by a site index x with sites within the interval 
[1, x] ([x + l,L]) belonging to region A (B) [thus we can 
label Sa by its site index, S(x)]. We stress that the QMC 
and DMRG results are on the same geometry and Hamil- 
tonian, and reproduce the same ground state energies; 
Figure 1 (and subsequent figures) should be considered 
as exact comparisons between S VB and 5 vN . 

The ID Heisenberg model is known to be critical and 
thus can be mapped to a 2D classical Hamiltonian at 
its critical point, which can be described by CFT in the 
limit L — > co. To address finite- length chains one can 
use the conformal mappings x — > x' — (L/ir) sminx/L) 
for PBC, and x — ► 2x' for OBC. Calculations within 
CFT [ll| obey S* vN (x) = (c/3) ln(a;') + Si for PBC, and 
5 vN (x) = (c/6)ln(2x') + ln( 5 ) + 5i/2 for OBC, where c is 
the central charge of the CFT, Si is a mo del- dependent 
constant, and g is Ludwig and Affleck's universal bound- 
ary term (l2T |. 
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FIG. 1: (color online) Entanglement entropies for a ID 
Heisenberg chain with PBC and OBC. Upper panels show 
the entropies as a function of the conformal distance x' = 
(L/n) sin(irx/L) for 100-site chains. Lower plots show the 
central charge c, obtained by fitting the numerical data to 
the CFT result, for several L. For PBC, c is calculated with 
the two smallest x' points removed. For OBC, the fits depend 
on the number of sites included, z, which we systematically 
decrease by removing x data points from the outside ends of 
the chain, c is shown for S' vN (closed symbols) and S VB (open 
symbols) for system sizes L = 64 (circles), L = 100 (squares), 
L — 128 (diamonds), and L — 200 (triangles) 



Figure Q] illustrates simulation results in both cases, 
the left (right) panels corresponding to PBC (OBC). For 
PBC both S YB and S vN appear to fit well to the CFT 
result, although S vN > S VB . The regression fit for S vN 
shows very good convergence with the central charge pre- 
dicted by CFT, c — 1, while the fit for S VB yields a lower 
value of c than predicted. For OBC both S vN and S VB 
split into two branches, the upper (lower) correspond- 
ing to an odd (even) number of lattice sites in A. This 
reflects a well-known "dimerization" effect induced by 
OBC Notice that contrary to the PBC case, now 

5 vN < S yB . A regression fit of the lower branch to the 
form (c/6) \n(2x') (inset) shows excellent convergence of 
S vN to the central charge predicted by CFT, c = 1, once 
finite-size effects and the proximity of the data to the 
open boundaries are taken into account. In contrast, 5 VB 
deviates significantly from the CFT result, giving c > 1 
when all or most data is included in the fit, changing to 
c < 1 as data closest to the open boundary is systemati- 
cally excluded [3], e.g. for z — L/2, cl-»oc ~ 0.85. 

Multi-leg ladders - Moving away from the ID chain, 
one can add "legs" to the lattice in a systematic way. In 
this case the sum over nearest neighbors is extended to 
neighbors along rungs as well as along legs. As noted 
before, DMRG imposes constraints to the subregion ge- 
ometry. In multi-leg ladders we choose to sweep in a ID 
path that visits first bonds sitting in rungs rather than 
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FIG. 2: (color online) Entanglement entropies for 3-leg (left) 
and 4-leg (right) ladder systems with OBC and 100 sites per 
leg. For odd-leg ladders, S(x) oc ln(x'). The left inset shows 
S(x) as a function of the conformal distance, x' , on a log 
scale. For even-leg ladders, S(x > £) = const. The right inset 
shows the site indexing used for multi-leg ladders where the 
bipartition A is shaded and labeled by x — 9. 




FIG. 3: (color online) Entanglement entropies divided by N, 
for JV-leg Heisenberg (filled symbols) and free-fermion (open 
diamonds) ladders, taken such that the region A includes 2N 2 
sites. For the Heisenberg model and large N, S VB oc NlnN, 
whereas S vN oc N. Data for free fermions, Sg N oc A" In N, are 
shown for comparison. We show data for ladders with length 
(sites per leg) L = 100 and L = AN. 



bonds sitting in legs (see Fig. [5]). DMRG computational 
demands increase dramatically with the number of legs, 
so in this paper we restrict ourselves to ladders with OBC 
up to N — 7 legs. QMC lacks this limitation and one can 
go up to N — 20 with minimal CPU effort. 

Figure [2] shows S vN and 5 VB calculated for the 3-leg 
and 4-leg ladder. As for the OBC ID chain, 5 VB > S vN . 
Entropy shows different behavior depending on N being 
even or odd. Even-leg ladders have a spin-gap [FoT ]. and 
thus only sites within distances from the boundary be- 
tween A and B smaller than the correlation length £ con- 
tribute to the entanglement, yielding S(x > £) = const. 
In contrast, odd-leg ladders are gapless, and thus all sites 
contribute to the entanglement, yielding S(x) oc ln(x'), 
which follows the CFT result in analogy to the ID case. 
As can be seen, S(x) splits into branches, with a (quasi-) 
periodic structure superimposed over the main depen- 
dence on x, the period being N (2N) for even- (odd-) leg 
ladders. This reflects the periodicity of the underlying ID 
path through which the algorithm sweeps, and the fact 
that valence bonds within the same rung are energetically 
favored 15j. The doubled period for odd-leg ladders is 
due to the same dimerization effect as in Ref. (l3j |. 

Area law in multi-leg ladders - We can use these results 
to address the question of the adherence of the 2D Neel 
state to an area law. To do so, we define the lattice geom- 
etry such that region A is rectangular, cutting a multi-leg 
ladder cleanly across all legs, so that the "area" separat- 
ing region A and B is equivalent to the number of legs 
in the ladder N. We choose the region A to contain 2N 2 
sites, to have an aspect ratio of order unity. In contrast, 
the entanglement entropy of a long narrow region would 
be dominated by the behavior of the gapless mode for 
odd- leg ladders. The 2:1 aspect ratio makes the region 



length even for all N, reducing even-odd oscillations. 

Figure[3]illustrates the simulation results for AT-leg lad- 
ders. Plotting S(x) /N versus N on a log scale, we obtain 
a multiplicative logarithmic correction to the «S VB area 
law, in agreement with results from Refs. 0, S|- How- 
ever, the linear slope is not present in the plot of the S vN 
data from the DMRG, which convincingly approaches a 
constant for large N [l6|. Clearly, for 5 vN the area law 
is indeed obeyed in the Neel groundstate, leading one 
to conclude that the multiplicative logarithmic correc- 
tion occurs in S' VB only. One can compare S VB to data 
obtained for free fermions, which has a well-known [ItJ 
logarithmic correction to the area law for S vN (Fig. [3]). 
In the next section, we explain S VB in the context of the 
bond length distribution in the QMC. We note also that, 
contrary to the suggestion in Ref. 0] , a gapless Goldstone 



mode will not give a logarithmic divergence to 5 vN since 



a gapless bosonic mode in 2D obeys an area law [18( . We 
have also done a spin-wave calculation of 5 vN for this sys- 
tem and found an area law, albeit with S(x)/N « 0.2, 
slightly lower than suggested for spin 1/2 in Fig. [3] 

Bond Length Distribution.- Sandvik defined the bond 
length distribution P(x,y) as the probability of a bond 
going from site x to site y, and found that P{x,y) ~ 
\x — y\~ p with p fts 3 in the Neel state [§]. This value of p 
gives the logarithmic divergence in as can be found 
by directly calculating 



cVB 
<~>4 



£ P(x,y)\n(2). 

x£A,y£B 



(3) 



We can understand the value of p from a scaling ar- 
gument similar to that in Ref. Consider the num- 
ber of bonds of length I exiting a region of linear size I. 
For p > 4, this scales to zero and such a state should 
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have no long-range order. For p = 4, this number is 
^-independent, corresponding to a critical state (p = 2 
is the critical power in ID, matching the observed loga- 
rithmic behavior of S VB ). For p < 2, the bond length 
distribution is unnormalized, and all bonds are long, so 
that there is no short-range order, p = 3 represents a 
state with both long- and short-range order. There is 
a heuristic argument that p — 3 corresponds to a state 
with Goldstone modes: the correlation function of two 
spins is the probability that they lie on the same loop. 
This can be estimated by the Green's function of a Levy 
flight with power law distributed steps with power 3, re- 
producing the power law decay of correlations. 

Discussion - In this paper, we have compared scal- 
ing properties of the valence bond entanglement entropy 
(S VB ) 0, [1| to the von Neumann entanglement entropy 
(5 vN ) in the spin 1/2 Heisenberg model on multi-leg lad- 
der geometries, using QMC and DMRG simulations. In 
ID, we find that S VB mimics the behavior of 5' vN closely, 
although it is less than 5 vN for periodic chains, and 
greater than S vN for open chains. In addition, fits to 
ID conformal field theory, which are excellent for S vN 
calculated via DMRG, appear to deviate significantly for 
S' VB in the large chain-size limit, approaching c < 1 for 
both boundary conditions 14 1. 

The fact that S* VB can be either greater or less than 
5 vN can be understood through simple examples. Let 
\(ij)(kl)...) denote a state in which sites i,j are in a 
singlet, sites k, I are in a singlet, and so on. Consider 
an 8-site chain, with sites 1 to 4 in region A. Then, 
the state | (12) (34) (56) (78)) + 1 (14) (32) (58) (76)) has van- 
ishing S\ B since no bonds connect A to B, but non- 
vanishing w 0.325. On the other hand, consider 
a 4-site chain, with sites 1 and 3 in region A. Then, the 
state |(12)(34)) + |(14) (32)> has a maximal S\ B , equal 
to 21n(2) w 1.386, while = ln(3) « 1.099 is smaller. 
This second state has the maximum possible Neel order 
parameter: it is the equal amplitude superposition of all 
configurations of bonds connecting the two sublattices. 
Thus, it is unsurprising that states with Neel order show 
5 VB > S" vN , a fact which we have demonstrated numeri- 
cally on multi-leg ladder systems with open boundaries. 

Defining the boundary between the two entangled re- 
gions as being bipartitioned by a cut across all legs on 
a ladder, we have shown using DMRG that 5 vN obeys 
the area law in the many-leg limit. Since DMRG can 
also accurately measure logarithmic size-dependencies of 
5 vN (in ID critical systems), this suggests that simula- 
tion procedures similar to those here might enable the 
measurement of area-law corrections in S* vN as indica- 
tors of exotic phases in other models, such as those with 
a spinon Fermi surface [ioj. 

The valence bond entanglement entropy harbors a mul- 
tiplicative logarithmic correction for the Neel ground- 
state, which we have shown is caused by the valence- 
bond length distribution, and is not present in 5* vN . It 



is clear that 5 VB is a reasonable measurement of entan- 
glement, readily accessible to numerical simulations in 
2D and higher, and capable of reproducing the area law 
in some gapped groundstates @, @]. However, the in- 
ability of S VB to provide a bound on S vN (unlike other 
measures such as Renyi entropies), along with its discrep- 
ancies from 5 vN in ID critical systems and the 2D Neel 
state, must be taken into account in proposals to use S' VB 
for future tasks such as characterizing topological phases 
or studying universality at quantum phase transitions. 
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